x
1
begin2
anim = for i in range(10, n_steps, step=10)3
i = min(i, n_steps)4
generate_plot(sol, (vals), 1:i)5
end6
gif(anim, "anim_fps15.gif", fps = 10)7
endx
1
begin2
anim2 = for i in range(10, n_steps, step=10)3
i = min(i, n_steps)4
generate_plot(sol2, (vals2), 1:i, τ=-τ)5
end6
gif(anim2, "anim2_fps15.gif", fps = 10)7
endx
1
generate_plot(sol, (vals), :)x
1
generate_plot(sol2, (vals2), :, τ=-τ)generate_plot (generic function with 1 method)x
1
function generate_plot(sol, vals, indices; τ=τ)2
(3
mean_values, below_mean, rolling_mean, 4
mean_rescaled_values, max_rescaled_values, min_rescaled_values,5
median_rescaled_values, q90_rescaled_values, q10_rescaled_values,6
sol_rescaled_matrix, rand_trajectory_indices7
) = vals8
l = [a; b; c]9
#indices = 1:i10
x = sol.t[indices]11
#p1 = plot(sol, yscale= :log, color=:gray, alpha=0.3, label="")12
p1 = plot(x, mean_values[1,indices], color=:red, label="mean")13
plot!(p1, legend=:outertopright, title="Mean wealth (τ=$τ, n=$n)")14
15
16
p3 = plot(x, below_mean[1,indices], label="below_mean")17
plot!(p3, x, rolling_mean[indices], label="time_avg")18
plot!(p3, legend=:outertopright, title="% below mean value")19
20
21
#p2 = plot(sol.t, sol_rescaled_matrix, color=:gray, alpha=0.3, label="")22
23
p2 = plot(x, mean_rescaled_values[1,indices], color=:red, label="mean")24
plot!(p2, x, max_rescaled_values[1,indices], color=:blue, label="max")25
plot!(p2, x, min_rescaled_values[1,indices], color=:blue, label="min")26
plot!(p2, x, median_rescaled_values[1,indices], color=:blue, label="50%")27
plot!(p2, x, q90_rescaled_values[1,indices], color=:magenta, label="90%")28
plot!(p2, x, q10_rescaled_values[1,indices], color=:magenta, label="10%")29
for i in rand_trajectory_indices30
plot!(p2, x, sol_rescaled_matrix[i,indices], color=:grey, label="")31
end32
33
plot!(p2, legend=:outertopright, title="Rescaled wealth") 34
plot(p1, p2, p3, size=(1000, 600), layout = l)35
end1×1002 Array{Float64,2}:
1.0 1.00441 1.00487 1.00757 1.00717 … -56820.3 -69083.6 -54654.5 -54654.51×1002 Array{Float64,2}:
0.0 0.506 0.527 0.514 0.516 0.526 0.515 … 0.419 0.418 0.414 0.417 0.4171×1002 Array{Float64,2}:
1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 … 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.01×1002 Array{Float64,2}:
1.0 1.12501 1.18447 1.25047 1.33236 … 7586.44 6779.1 8662.02 8662.021×1002 Array{Float64,2}:
1.0 0.828039 0.798682 0.786456 … -1585.44 -1305.03 -1625.82 -1625.821×1002 Array{Float64,2}:
1.0 0.999422 0.997064 0.996852 0.99739 … -3.2702 -2.77611 -3.59646 -3.596471×1002 Array{Float64,2}:
1.0 1.05807 1.08189 1.10222 1.11794 … 74.9351 63.5709 82.4535 82.45361×1002 Array{Float64,2}:
1.0 0.943985 0.92129 0.907465 0.8931 … -89.409 -75.4356 -96.7221 -96.72221000×1002 Array{Float64,2}:
1.0 0.906762 0.916616 0.898629 0.884472 … 1.24771 1.59119 1.59119
1.0 0.945487 0.858075 0.862511 0.784747 101.542 135.162 135.163
1.0 0.941019 0.873183 0.851367 0.803476 14.8362 18.6452 18.6452
1.0 0.946316 0.928038 0.905564 0.883512 14.0832 17.7311 17.7312
1.0 0.97697 0.921472 0.832453 0.915906 36.9134 46.7759 46.776
1.0 0.987504 0.889869 0.896703 0.892412 … 15.576 19.3141 19.3141
1.0 0.999409 1.0398 1.04312 1.10827 -30.0248 -37.1702 -37.1702
⋮ ⋱ ⋮
1.0 0.946462 0.958537 0.961702 0.960493 51.1245 66.7855 66.7855
1.0 0.981512 0.979709 0.907478 0.91932 … 453.653 599.627 599.627
1.0 0.936449 0.95218 0.968702 0.945949 -1.0746 -1.40474 -1.40474
1.0 1.03358 1.04421 1.06579 1.06526 -8.57554 -11.2261 -11.2261
1.0 1.03495 1.05807 1.01143 0.990378 -18.5972 -22.4368 -22.4368
1.0 1.00405 1.00977 1.02101 0.992355 -21.5162 -29.1503 -29.1503x
1
1
n_steps, vals = get_vals(sol)1×1002 Array{Float64,2}:
1.0 1.00028 1.00226 1.00322 1.00386 … 7.22464 7.25642 7.28517 7.285171×1002 Array{Float64,2}:
0.0 0.497 0.52 0.523 0.516 0.516 0.531 … 0.565 0.568 0.562 0.572 0.5721×1002 Array{Float64,2}:
1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 … 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.01×1002 Array{Float64,2}:
1.0 1.14538 1.22133 1.28444 1.30007 … 2.33007 2.3128 2.43117 2.431171×1002 Array{Float64,2}:
1.0 0.858851 0.807647 0.779717 … 0.45922 0.465704 0.464225 0.4642251×1002 Array{Float64,2}:
1.0 1.00044 0.997208 0.996108 … 0.958374 0.963146 0.963132 0.9631321×1002 Array{Float64,2}:
1.0 1.05764 1.08443 1.0969 1.11043 … 1.34074 1.33863 1.33335 1.333351×1002 Array{Float64,2}:
1.0 0.942907 0.922381 0.90541 … 0.713265 0.715867 0.713505 0.7135051000×1002 Array{Float64,2}:
1.0 0.997736 0.98027 0.982883 … 0.873308 0.848167 0.791675 0.791676
1.0 0.950797 0.955326 0.968567 0.722607 0.723182 0.713573 0.713573
1.0 1.01378 0.990531 0.963829 0.729599 0.733584 0.712423 0.712424
1.0 0.95753 0.936931 0.961542 0.737401 0.713991 0.683079 0.683079
1.0 0.918856 0.965992 0.963822 1.51141 1.57924 1.61247 1.61247
1.0 0.955954 0.996973 1.00992 … 1.09496 1.09423 1.07979 1.07979
1.0 0.964941 0.973986 0.985951 0.768035 0.758382 0.770829 0.770829
⋮ ⋱ ⋮
1.0 1.02924 0.956967 0.955749 1.14021 1.17685 1.18437 1.18437
1.0 1.03407 1.045 1.05124 … 0.996936 0.963976 0.999022 0.999021
1.0 1.08054 1.11138 1.13716 0.657496 0.652194 0.653937 0.653936
1.0 1.03727 1.01601 1.03225 1.3156 1.30934 1.31979 1.31979
1.0 1.04184 1.04809 1.03317 0.99756 1.02491 1.11198 1.11198
1.0 1.05685 1.09597 1.11517 0.567786 0.605907 0.581047 0.581047x
1
n_steps2, vals2 = get_vals(sol2)get_vals (generic function with 1 method)x
1
function get_vals(sol)2
sol_matrix = hcat(sol.u...)3
mean_values = mean(sol_matrix, dims=(1))4
n_steps = length(sol.t)5
below_mean = mean(sol_matrix .< mean_values, dims=1)6
rolling_mean = cumsum(below_mean, dims=(2))[1,:] ./collect(range(1.0, n_steps, step=1.0))7
sol_rescaled = deepcopy(sol)8
sol_rescaled_matrix = sol_matrix ./ mean_values9
#sol_rescaled.u = [r[:] for r in eachrow(sol_rescaled_matrix)]10
mean_rescaled_values = mean(sol_rescaled_matrix, dims=(1))11
max_rescaled_values = maximum(sol_rescaled_matrix, dims=(1))12
min_rescaled_values = minimum(sol_rescaled_matrix, dims=(1))13
median_rescaled_values = median(sol_rescaled_matrix, dims=(1))14
q90_rescaled_values = hcat([quantile(r[:], 0.9) for r in eachcol(sol_rescaled_matrix)]...)15
q10_rescaled_values = hcat([quantile(r[:], 0.1) for r in eachcol(sol_rescaled_matrix)]...)16
rand_trajectory_indices = rand(1:n_steps, 5)17
vals = (18
mean_values, below_mean, rolling_mean, 19
mean_rescaled_values, max_rescaled_values, min_rescaled_values,20
median_rescaled_values, q90_rescaled_values, q10_rescaled_values,21
sol_rescaled_matrix, rand_trajectory_indices22
)23
return n_steps, vals24
end[36mSDEProblem[0m with uType [36mArray{Float64,1}[0m and tType [36mFloat64[0m. In-place: [36mfalse[0m
timespan: (0.0, 100.0)
u0: [1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0 … 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0]x
1
begin2
μ=0.0213
σ=0.144
τ=-0.155
n=10006
tspan=(0.0, 100.0)7
prob = get_diff_eq(μ=μ, σ=σ, τ=τ, n=n, tspan=tspan)8
end[36mSDEProblem[0m with uType [36mArray{Float64,1}[0m and tType [36mFloat64[0m. In-place: [36mfalse[0m
timespan: (0.0, 100.0)
u0: [1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0 … 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0]x
1
prob2 = get_diff_eq(μ=μ, σ=σ, τ=-τ, n=n, tspan=tspan)xxxxxxxxxx1
1
using DifferentialEquationsxxxxxxxxxx1
1
using Plots; pyplot() #plotly() # Using the Plotly backendxxxxxxxxxx1
1
using StatisticsFloat64xxxxxxxxxx1
1
typeof(1/2^(4))get_diff_eq (generic function with 1 method)xxxxxxxxxx14
1
function get_diff_eq(;2
μ::Number=1.0, 3
σ::Number=1.0, 4
τ::Number=0.0, 5
n::Number=2, 6
tspan::Tuple{Number,Number}=(0.0, 100.0), 7
#u₀::Array{Number}=nothing8
)9
u₀ = ones(n)10
f(u,p,t) = μ*u - τ*(u .- mean(u))11
g(u,p,t) = σ*u12
prob = SDEProblem(f,g,u₀,tspan)13
return prob14
endretcode: Success
Interpolation: 1st order linear
t: 1002-element Array{Float64,1}:
0.0
0.1
0.2
0.30000000000000004
0.4
0.5
0.6
⋮
99.59999999999862
99.69999999999861
99.7999999999986
99.8999999999986
99.9999999999986
100.0
u: 1002-element Array{Array{Float64,1},1}:
[1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0 … 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0]
[0.998019015752889, 0.9510658414926022, 1.0140694777315076, 0.9578012096701445, 0.919115750512006, 0.9562252003478064, 0.9652141479108769, 0.9706803007418675, 1.0693271786638994, 1.005725074125227 … 1.0260034142387402, 0.9774434436720234, 1.0422001114652972, 1.0055480480456966, 1.0295297866000548, 1.0343606289661411, 1.0808454876960563, 1.037564871012976, 1.0421317418173428, 1.0571462217470442]
[0.9824844006466523, 0.9574841946060301, 0.9927687303308445, 0.9390475395041075, 0.9681739420261597, 0.9992250635127821, 0.9761864988499609, 0.9749777440755799, 1.1185920911819842, 0.9996010984096022 … 1.0128148913950608, 0.9643601311596616, 1.0871619763929496, 1.0071108110479834, 0.9591293404348149, 1.0473567576765512, 1.1138866594668928, 1.0183030903754369, 1.0504537398543374, 1.098447284362764]
[0.9860463848081689, 0.9716847501087723, 0.9669308623772892, 0.9646365640051681, 0.9669246811117161, 1.0131717110854639, 0.9891246079945577, 1.0543268531059922, 1.1732885626077654, 1.0152492187419375 … 1.0129504149464972, 0.9727387036196269, 1.0630960079494067, 1.0176841188785029, 0.9588252695398269, 1.0546222263836145, 1.140818110044002, 1.0355764384147492, 1.036495271292124, 1.1187574179226343]
[1.0085711020573354, 0.9688124576973165, 0.9731906203406078, 0.8875361566189837, 0.9337511221653139, 0.9917939300864601, 0.9711477516182847, 1.1034427215447837, 1.2381125012644685, 0.9718812833537678 … 1.053858781600805, 1.0088169286450392, 1.0233671401856204, 1.0105631872778331, 0.986423029088711, 1.106337284755494, 1.0571880222513896, 1.0279559033203622, 1.062273669287912, 1.0744752535914042]
[0.9985448018002562, 0.9890005186576608, 0.9814689545753578, 0.8657904412760288, 0.9652600449580027, 0.9654300719754308, 0.9410039754317325, 1.0909291826622525, 1.2736813425671867, 1.0005982339652475 … 1.0331588937786556, 0.9534195061940692, 1.0505014177103604, 1.0226781364606887, 0.9381766550748648, 1.0908547230437289, 0.9666615732892213, 1.0158855551936894, 1.1867798328282506, 1.141474031208177]
[0.9979176373268588, 0.9677737210159575, 0.9731578453695062, 0.8386383469588713, 0.9092004841646534, 0.9653225608363216, 0.9004746679142716, 1.116162018784377, 1.226907049852087, 0.9847765727788872 … 0.994266256117669, 0.9831930560695477, 1.0838941959280866, 0.985100017865773, 0.9028551378884836, 1.1072539076277403, 0.9697196118730311, 1.0838009688192383, 1.1887551374310195, 1.1301885781373155]
⋮
[6.007617830026853, 5.518750059680582, 5.208741916089139, 4.9526517022883745, 10.717076559978665, 7.178555965772261, 5.78486641092313, 5.846649048401186, 7.127984106953892, 7.74830558098805 … 4.628301277046172, 4.9016179248704, 7.077297282476337, 6.589564669243558, 8.53420124280035, 6.671835660288932, 4.94988196662736, 8.904778362042173, 7.323413010417348, 4.287439123817096]
[6.173959160967872, 5.364251369804051, 5.247376240590091, 5.099542032946705, 11.094601933625011, 7.764965030854282, 5.529084425099932, 6.040577700122486, 7.250356506418934, 7.811905893952469 … 4.863641853013579, 4.774401730694823, 6.937836279355249, 6.853764932166131, 7.995241148298569, 7.0448134838736935, 4.8801853329138964, 9.157465845073196, 7.06475436088209, 4.369006478703452]
[6.30933043468796, 5.220569794041702, 5.2710862516930135, 5.327456598058746, 10.919420068282582, 7.910723182676571, 5.548771339988939, 6.339042836253449, 6.720111885860936, 7.5351288104820675 … 4.594848297650571, 4.691698618421832, 6.940738679965065, 7.0549749306087195, 8.237599513823463, 7.202496228653328, 4.750167242943078, 9.504738438020553, 7.207004525988839, 4.102050277005823]
[6.154656444749044, 5.247713695808727, 5.323191087075731, 5.181014277902353, 11.459615652070859, 7.940175370971785, 5.503139422508213, 6.797889948885585, 6.8397785043376915, 7.659744915814635 … 4.68146931874797, 4.7683275457256, 7.416770606550524, 6.847848253863965, 8.539712640187402, 6.995013967970305, 4.732595626015025, 9.50109794046614, 7.437169824974048, 4.396714226900677]
[5.767490375177858, 5.198504465040368, 5.190126802966306, 4.976346772957234, 11.74710704521954, 7.866463333387993, 5.615621220436422, 6.6074258655017175, 6.516238283967131, 8.076917370797762 … 4.9988256876993535, 4.96913448861327, 7.18312593044942, 6.9249006705389204, 8.628344057868167, 7.278043384656616, 4.764039686985904, 9.61492984870188, 8.101000602930913, 4.233026003385014]
[5.7674920298146075, 5.198504302704035, 5.19012800986175, 4.976346616834273, 11.747106687349568, 7.86646181150715, 5.615620287656042, 6.607425348017516, 6.516237627450753, 8.076914895411957 … 4.998824863734882, 4.969134897876926, 7.183125910726612, 6.924900610903935, 8.628345711685018, 7.278041321661102, 4.764038905983389, 9.614930051574952, 8.10100132320552, 4.233025305700403]x
1
begin2
dt=0.13
sol = solve(prob,EM(),dt=dt)4
sol2 = solve(prob2,EM(),dt=dt)5
endArray{Float64,2}xxxxxxxxxx1
1
typeof(hcat(sol.u...))xxxxxxxxxx1
1
rand(1:n_steps, 10)xxxxxxxxxx6
1
html"""<style>2
main {3
margin: auto;4
max-width: 1400px;5
}6
"""